Quasi-stationary simulations of the contact process on quenched networks 
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We present high-accuracy quasi-stationary (QS) simulations of the contact process in quenched 
networks, built using the configuration model with both structural and natural cutoffs. The critical 
behavior is analyzed in the framework of the anomalous finite size scaling which was recently shown 
to hold for the contact process on annealed networks, ft turns out that the quenched topology does 
not qualitatively change the critical behavior, leading only (as expected) to a shift of the transition 
point. The anomalous finite size scaling holds with exactly the same exponents of the annealed 
case, so that we can conclude that heterogeneous mean-field theory works for the contact process 
on quenched networks, at odds with previous claims. Interestingly, topological correlations induced 
by the presence of the natural cutoff do not alter the picture. 
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I. INTRODUCTION 

Dynamical processes are strongly affected by the struc- 
ture of the pattern which mediates interactions. This 
common wisdom, whose evidence has been confirmed 
by decades of investigations in statistical mechanics, has 
assumed a new relevance in recent years, with the ad- 
vent of complex networks theory. Percolation, epidemic 
spreading or synchronization are among the phenomena 
for which the effect of a disordered topology has been 
more thoroughly studied [IH3], but for many other types 
of dynamics a complex substrate structure also induces 
novel and nontrivial dynamical behavior [4j [5] . This oc- 
curs in particular when the pattern of connections is ex- 
tremely heterogeneous, as in scale-free (SF) networks, 
where the probability that an element is connected to k 
others (the degree distribution) is given by a power-law 
form P(k) ~ A: -7 0. 

Along with the recognition that dynamics on networks 
can be very different from lattices, comes the natural 
question about how theoretical methods can be adapted 
to deal with topologically complex substrates. It has 
been early recognized that the standard mean-field ap- 
proach must be modified on strongly heterogeneous net- 
works, to take into account the broad variability in the 
connectivity of vertices [2J. Thus, in the heterogeneous 
mean-field (HMF) theory [4j [5] the order parameter is 
replaced by a set of analogous quantities which depend 
explicitly on the degree k of the node considered (the de- 
gree of a node is the number of other nodes to which 
it is directly linked). This simple modification turns 
out to describe with remarkable accuracy the behavior 
of many systems, reinforcing the naive expectation that 
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mean-field methods must work on networks, due to their 
infinite-dimensional nature. 

For this reason it came as a surprise when the numer- 
ical investigation of the behavior of the contact process 
(CP) on generic power-law distributed networks gave re- 
sults in apparent disagreement with the predictions of 
HMF theory [7]. The contact process [8] is an extremely 
simple model for spreading, which, in a general network, 
is defined as follows: Vertices can be in two different 
states, either empty or occupied. The dynamics includes 
the spontaneous annihilation of occupied vertices, which 
become empty at unitary rate, and the occupation of 
empty neighbors by occupied vertices, with rate A/fcj, 
where ki is the degree of the occupied node. The model 
is characterized by a phase transition at a value of the 
control parameter A = A c , separating an active phase 
from an absorbing one, devoid of occupied vertices [5]. 

In Ref. [7] the HMF theory for CP was derived in the 
limit of infinite network size. Its predictions could not be 
directly checked against numerical simulations because 
of the presence of extremely large finite size effects. A 
comparison was made possible by the introduction of a 
finite size scaling (FSS) ansatz [TU], adapted to a network 
topology, leading to the conclusion that CP dynamics 
on quenched networks was not described by the HMF 
approximation. 

Such a claim was criticized by Park and collabora- 
tors [TTJ [12], which later proposed an alternative FSS 
ansatz, based on a droplet excitation theory [T3]. A dis- 
crimination between the two approaches turned out to be 
nontrivial: Even for simulations on annealed networks, 
which are expected to be described exactly by mean-field 
theory, numerical results did not satisfactorily conform 
to any of the two competing theories 12,. 

Clarifying results on this issue have been obtained 
later [TJ1 \T5\ , revealing that FSS on annealed networks 
is actually more complicated than previously assumed. 
At odds with what happens in lattices, the behavior of 
the CP on networks of finite size depends not only on 



the number of vertices N but also on the moments of 
the degree distribution. This last feature implies, for 
SF networks, that the scaling around the transition de- 
pends explicitly on how the largest degree k c diverges 
with N. Such a dependence (which both previous FSS 
approaches were lacking) introduces very strong correc- 
tions to scaling. However, if such corrections are properly 
taken into account, it is possible to show that the CP on 
annealed networks agrees, with high accuracy, with the 
predictions of HMF theory [16 . A detailed comparison 
between HMF and precise numerical results for annealed 
networks has been possible by using the quasi-stationary 
(QS) state simulation method [TTJ, [IS]. This is an op- 
timized numerical technique which prevents the system 
from falling into the absorbing configuration, thus allow- 
ing the detailed investigation of the order parameter close 
to the transition for finite systems. 

In the present paper, we go beyond the results of 
Ref. [16] and exploit QS simulations for studying the CP 
transition on quenched SF networks. We assume that 
finite size scaling on quenched networks has exactly the 
same anomalous form valid on annealed ones. In this way 
we are able to determine, with an unprecedented accu- 
racy, the location of the critical point and the value of 
the exponents associated with the transition. It turns out 
that the exponents computed numerically for quenched 
networks are in excellent agreement with the HMF pre- 
dictions: The CP on networks definitely obeys heteroge- 
neous mean-field theory. This result fills the gap in the 
debate about the validity of HMF theory for CP. Both 
competing theories [7J [13] were not correct. The appar- 
ent discrepancy between theory and numerics reported in 
Ref. [7] was due to an incorrect FSS theory and the pres- 
ence of large corrections to scaling in simulations. The 
alternative FSS approach of Ref. [13] was incorrect as 
well. The correct FSS theory is more complicated than 
previously assumed, but perfectly accounts for numerical 
results on both annealed and quenched networks. 

The present paper is structured as follows: In Sec- 
tion II we describe the quenched networks considered, 
the contact process dynamics, and the QS simulation 
method. The numerical determination of the transition 
point and the associated exponents are presented in Sec. 
III. Our concluding remarks are discussed in Sec. IV. 



II. MODEL AND SIMULATION METHODS 

The CP is simulated on quenched networks with N ver- 
tices and degree sequence {fei, k%, . . . ,fcj\r}. The networks 
are built according to the configuration model (CM) [T§] . 
where the degree of each vertex is defined at the begin- 
ning of the simulation as a random variable with a power 
law distribution P(k) = Ak 1 and ko < k < k c . The 
degree distribution has a hard degree cutoff k c — N 1 ^ 
(degrees larger than k c are not allowed) , where uj > 7 — 1 
is the cutoff exponent. This choice reduces sample to 
sample fluctuations that make the critical analysis dif- 



ficult even for the simpler case of annealed networks 

palm eq]. 

Two models for the link assignment are investigated, 
namely, random and ordered configuration models. In 
the random CM [21] , two stubs (not connected links) are 
selected at random and the respective nodes connected, 
avoiding self and multiple connections. The procedure is 
repeated until all stubs are joined. In the ordered CM 
nodes are sorted in descending degree order. Then a 
stub of the most connected node is connected to another 
randomly selected stub, avoiding self and multiple con- 
nections. After all stubs of the largest vertex are con- 
nected, the next most connected node with free stubs is 
picked up and the procedure is repeated until all stubs are 
connected. The ordering of nodes allows to build large 
networks even for 7 — 1 < u) < 2, which is usually impos- 
sible with random selection of nodes. A consequence of 
this nonrandom procedure is that disassortative degree 
correlations are generated [2"2"] . 

The CP simulations on an arbitrary network are per- 
formed with the standard protocol [5]: At each time step, 
an occupied vertex j is chosen at random and time is up- 
dated as t ->• t + At, where At = 1/[(1 + X)n(t)} and n{t) 
is the number of occupied vertices at time t. With proba- 
bility p = 1/(1 + A), the occupied vertex becomes vacant. 
With complementary probability 1 — p = A/(l + A), one 
of its kj neighbors is randomly selected and, if empty, 
it becomes occupied. If the selected neighbor is already 
occupied nothing happens and the simulation proceeds 
to the next step. 

The standard numerical procedure to investigate the 
finite size scaling at absorbing phase transitions is based 
on the determination of the average of the order param- 
eter (in this case the density of active nodes), p s , re- 
stricted only to surviving runs. Such a technique is quite 
inefficient, because surviving configurations are very rare 
at long times, and in order to get precise results an ex- 
ceedingly large number of realization of the process is 
needed. An alternative strategy consists in constraining 
the system in a quasi-stationary state and measuring its 
properties. In practice this is implemented by replacing 
the absorbing state, every time the system tries to visit 
it, with an active configuration randomly taken from the 
history of the simulation [T7J . For this task, a list of M 
active configurations is stored and constantly updated. 
An update consists in randomly choosing a configuration 
in the list and replacing it by the present active con- 
figuration with a probability p r At. After a relaxation 
time t r , the QS quantities are determined during an av- 
eraging time t a . The QS probability P n that n vertices 
are occupied is computed during the averaging interval, 
each configuration with n active vertices contributing to 
the QS distribution with a probability proportional to its 
lifespan, i.e. proportional to 1/n. From the particle dis- 
tribution P n , the characteristic properties of the QS state 



can be computed, such as the density of active nodes 
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and the characteristic relaxation time [HjI [TH] 
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The QS method, which had been previously used to accu- 
rately determine the universality class of several models 
with absorbing configurations [17j [23] [24] , was recently 
applied to investigate the critical properties of the CP in 
annealed SF networks [16 . 

The FSS method for non-equilibrium absorbing phase 
transitions usually assumes, at the transition, an asymp- 
totic power-law dependence with the system size of the 
order parameter, p s ~ N~ v , and of the characteristic 
temporal scale, t ~ N a . Values of the control parame- 
ter above and below the critical point give rise to posi- 
tive and negative curvatures in plots of In p s versus In N, 
respectively, for large N [9]. The CP in annealed net- 
works has a different, anomalous, scaling form [T4l IT51 |2"U] 
which includes an explicit dependence on the factor g = 
(k 2 )/(k) 2 . In Ref. [TH] the same form was analytically 
shown to hold on annealed networks for the QS density 
p and characteristic time r, namely 
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The presence of the factor g, which can be determined 
exactly for a given degree distribution P(k), brings an 
additional dependence on N for 7 < 3. Hence, the expo- 
nents predicted by HMF theory for the critical properties 
as a function of N are HU : 



1 (Z-l n 

v = — h max , 

2 \ 2lj ' 

1 ( 3 -"< n 

a = 2" max l^T' 



III. RESULTS 
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The QS simulations are performed on networks with 
sizes ranging from N = 10 4 up to N = 2 x 10 7 and up 
to times t s = 2 x 10 6 . The so-called structural (ui = 2) 
and natural (ui = 7 — 1) hard cutoffs [25j are consid- 
ered. The structural cutoff leads to the absence of degree- 
degree correlations [25] in the random CM. In this case, 
the procedure to build the network corresponds to the 
uncorrelated configuration model (UCM) [5T]. On the 
contrary, the ordered CM has, for 7 < 3, disassortative 
degree-degree correlations, which are disregarded in the 
HMF analysis [U [5] . The determination of the QS dis- 
tributions starts after a relaxation time t r = 10 6 . The 




FIG. 1: (Color online) Moment ratios as function of the dis- 
tance to the critical point A = A — A c for CP on annealed SF 
networks with degree exponent 7 = 2.75 and cutoff exponent 
W = 2. Network sizes N = 8 x 10 4 , 16 x 10 4 , 32 x 10 4 , and 
64 x 10 4 are shown. 



dynamics is run for several independent network realiza- 
tions in order to average over many topological configu- 
rations. Between 200 and 500 network samples were used 
for each value of N. It is important to notice that the 
critical relaxation time is very short when compared to 
the same quantity for regular lattices below the critical 
dimension, due to the small world property of complex 
networks. 



A. Determination of the critical point 

In a quenched network the investigation of the transi- 
tion is made difficult by the presence of dynamical cor- 
relations, which shift the position of the critical point. 
From the analytical point of view the simplest way to 
take them into account is by means of a homogeneous 
pair approximation, which yields |26J A c = (k) /((k) — 1). 

The criterion of null curvature, usually applied to 
find numerically the critical point in regular lattices (see 
Sec.|H| was used for the CP on quenched networks [71 [T3"]. 
However, as shown for annealed networks [16J, such 
method for estimating the critical point can be mislead- 
ing if anomalous scaling of the form (13]) occurs. Addi- 
tionally, dynamical analyses like spreading experiments 
and density decay in time [9] require absence of finite size 
effects which are unavoidable in random networks with 
the small world property [IS] [27]. For these reasons, we 
determine the position of the critical point by analyz- 
ing the ratio between moments of the order parameter. 
The well-known fourth order Binder reduced cumulant 
Ui = 1— (p 4 )/3(p 2 ) 2 is a standard quantity to determine 
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FIG. 2: (Color online) Main plot: Moment ratios (p 2 )/(p) 2 
as a function of the annihilation probability for CP on UCM 
networks with degree exponent 7 = 2.75 and minimum de- 
gree fco = 6. Symbols represent stochastic simulations, while 
solid lines are hyperbolic tangent regressions, drawn as guides 
to the eyes. Inset: The corresponding fourth order moment 
ratios (p 4 )/(p 2 ) 2 around the critical point. Network sizes 
N = 1.28 x 10 6 , 2.56 x 10 6 , 5.12 x 10 6 , and 1.024 x 10 7 are 
shown. 



critical points in equilibrium phase transitions on mag- 
netic systems [25] • At criticality, this cumulant does not 
depend on size, implying that curves U4,(X;N) against 
the control parameter A (the temperature in magnetic 
systems) for different sizes all cross at A = A c . More in 
general, moment ratios defined as 



M n = JP± 



(6) 



are expected to be size-independent (and thus cross- 
ing each other) at the transition. Universal and size- 
independent moment ratios were studied for absorbing 
phase transitions in lattice models 29 31 . The size in- 
dependence of moment ratios in lattices systems results 
from the scaling invariance close to the critical point (see 
e.g., Ref. [32j for a recent review). The order parameter 
in the CP, the density of active vertices, is positive defi- 
nite and thus, in contrast with the magnetic counterpart 
for which all odd moments are null due to symmetry, 
there is no restriction on the values of q and s that can 
be considered. 

In order to probe the validity of this method for finding 
the critical point in complex networks, we first determine 
moment ratios for the CP in annealed SF networks, for 
which the critical point is exactly known to be A c = 1 
[15] . Figure[T|shows the moment ratios M" for CP on an- 
nealed networks with degree exponent 7 = 2.75 and dif- 
ferent sizes. As observed for macroscopic quantities like 
density or characteristic time in annealed networks |16j . 
moment ratios also have finite size corrections; the cross- 
ing points converge to A = 1 strictly only for N — > 00. 
The convergence is fast, particularly for the higher order 
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FIG. 3: (Color online) Third order moment ratios for CP 
on UCM networks with degree exponents: (a) 7 = 2.25, (b) 
7 = 2.50, and (c) 7 = 3. The minimum degree is fco = 
6. Symbols represent stochastic simulations, while solid lines 
are hyperbolic tangent regressions, drawn as guides to the 
eyes, (d) Extrapolation used to determine the critical point 



for 7 = 2.25. Network sizes TV 
1.024 x 10 7 are shown. 
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moment ratios, which, however, are more susceptible to 
large statistical fluctuations. For this reason we analyze 
moment ratios up to fourth order. We find that the cross- 
ing points of moment ratios tend to constant values, i.e. 
independent on the network degree distribution (data not 
shown). The second order moment ratio Mf-^ = 1.667(3) 
is slightly larger than the value 1.660 obtained for CP 
on the complete graph [IB] . The other moment ratios 
determined are M^ = 2.190(4) and Af 2 4 2 = 3.452(3), re- 
spectively. 

In the case of quenched networks, the moment ratios 
for CP on UCM networks (w = 2) with 7 = 2.75 and 
fco = 6 versus the annihilation probability p = 1/(1 + A) 
are plotted in Fig. [2] Moment ratios are determined 
for each network realization, and then averaged over the 
ensemble of networks considered. The simulation data 
are very accurate, allowing to clearly resolve the crossing 
points taking place between 0.47965 and 0.47970 for all 
moment ratios analyzed (including the third order which 
is not shown.) As a consequence the critical point posi- 
tion is estimated to be p c = 0.47968(3) or, equivalcntly, 
A c = 1.0847(1), where the number in parenthesis repre- 
sent the numerical uncertainty in the last digit. The crit- 
ical moment ratios are similarly estimated, yielding the 
values Mfj = 1.77(2), Mf x = 2.42(2) and M| 2 = 4.10(9). 

Third order moment ratios for different values of 7 
are shown in Fig. [3Fa)-(c). We observe that the higher 
the network heterogeneity, the stronger the finite size 
effects; thus larger systems are required for small 7 in 



Pc 



M?! 



M-; 



2.25 


0.4896(1) 0.48861 


1.92(1) 


2.77(3) 


5.1(1) 


2.50 


0.48425(5) 0.48386 


1.85(2) 


2.62(4) 


4.7(1) 


2.75 


0.47968(3) 0.47969 


1.77(2) 


2.42(4) 


4.10(9) 


3.00 


0.47602(3) 0.47628 


1.72(2) 


2.32(3) 


3.81(7) 


3.25 


0.47303(3) 0.47356 


1.70(1) 


2.25(2) 


3.60(7) 



ANN 1/2 



1/2 1.667(3) 2.190(4) 3.452(3) 



TABLE I: Critical points p c — 1/(1 + A c ) and critical mo- 
ment ratios for CP on UCM networks with minimum degree 
fco = 6. The third column reports the predictions of the pair 
approximation pj = ((fc) — l)/(2(fe) — 1). The values corre- 
sponding to annealed networks (ANN) are also included for 
sake of comparison. 



order to determine the critical point with little uncer- 
tainty. The critical points are determined by extrap- 
olating the crossing points between curves correspond- 
ing to sizes N a and Nf, against the inverse of the mean 
size 1/N = 2/(N a + Nt,), assuming a functional depen- 
dence \ C (N) = X c + const/TV, as shown in Fig. UMd) for 
7 = 2.25. The critical points and moment ratios com- 
uted for different degree exponents are shown in Table 
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It is worth to stress that if we replace N~ in the 
FSS by N~ b , where & is a fit parameter, the estimates 
vary inside estimated error bars. For the most extreme 
case, 7 = 2.25, we obtain b = 0.72 and the critical point 
changes to A c = 0.4897. Our simulations show that the 
moment ratio analysis is a very efficient method to de- 
termine the critical point of absorbing phase transitions 
in complex networks. The moment ratios, in opposition 
to annealed networks, depend now on the degree expo- 
nent through an exponential form: As 7 increases, they 
approach the value of the annealed case. The position 
of the critical point, on other hand, is remarkably close 
to the predictions of pair approximation calculations [26] 
and it approaches the annealed value p c = 1/2 for 7—^2 
when (k) — > 00. 



B. Critical exponents 

An accurate knowledge of the critical point allows to 
determine the exponents characterizing the singular be- 
havior at the transition and also to probe the presence 
of corrections to the scaling. Assuming that Eqs. ([3| 
hold also for quenched networks, we perform a power- 
law regression analysis in the plots of lnp vs IniV and 
Inp vs higN, as well as a similar analysis for the charac- 
teristic time r. Both are shown in Fig. |4J System sizes 
smaller than N = 10 5 are excluded from the regression to 
avoid subleading corrections to scaling. The exponents 
obtained numerically as well as the HMF predictions in 
Eqs. Q and (5) arc shown in Table [II] and Fig. [5] 

For 7 = 2.75, 3 and 3.25, the exponents and a esti- 
mated from a direct power law regressions as a function of 
N are quantitatively consistent with those measured for 
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FIG. 4: (Color online) Critical QS quantities (density in the 
top and characteristic time in the bottom) for the UCM net- 
work with minimum degree fco = 6 and distinct degree expo- 
nents. Symbols represent QS simulations and solid lines the 
respective power law regressions. Dashed lines have slopes 
±1/2. 



annealed networks [16] but they are larger and smaller, 
respectively, than the HMF values. This is due to the 
fact that HMF exponents [Eq. (pjj)] are calculated by con- 
sidering the asymptotic scaling of g, g ~ k^~^ . However, 
the most relevant corrections to scaling 



g ~ const x 
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7-2' 



fc c 3 ' 7 (7) 



are nonnegligible even for the largest system sizes consid- 
ered here. Notice that for both 7 — >■ 2 and 7 — > 3 correc- 
tions to scaling decay only logarithmically. If instead we 
perform the regressions p ~ (gN)~ Sv and r ~ {N/g) Sa 
we obtain exponents very close to 0.5, in full agreement 
with the HMF prediction Shmf = 1/2. 

For 7 = 2.50, the corrections to the scaling of factor 
g vanish as A^ -1 / 4 , the fastest decay in the analyzed 7 
range. Both characteristic time and density exponents 
agree with HMF theory if the factor g is explicitly in- 
cluded. If g is not included, the exponents still agree 
with the HFM predictions, due to the fast decay of the 



7 



VHMF 



OiHMF 



sv 



Sac Shu 



2.25 0.69(2) 0.6875 0.31(2) 0.3125 0.52(2) 0.43(3) 1/2 

2.50 0.63(1) 0.6250 0.37(1) 0.3750 0.51(1) 0.49(2) 1/2 

2.75 0.57(1) 0.5625 0.43(1) 0.4375 0.50(1) 0.50(2) 1/2 

3.00 0.55(1) 1/2 0.44(2) 1/2 0.51(1) 0.48(2) 1/2 

3.25 0.51(2) 1/2 0.47(2) 1/2 0.49(2) 0.49(2) 1/2 

TABLE II: FSS exponents for the CP on UCM network with fco = 6. HMF exponents, Eq. (pi), are also shown for comparison. 
The number in parenthesis represents the error in the last digit given by the standard deviation of the exponents fitted for the 
two curves that are closest to the critical point, being one below and other above A c . 
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FIG. 5: (Color online) FSS exponents for CP on UCM net- 
works with minimum degree fco = 6, compared with the HMF 
and droplet mean field (DMF) [13] theory predictions. Main 
panel shows the results of power law regressions p ~ N~ v 
and p ~ {gN)~ while the inset shows similar analysis for 
the characteristic time r. 




FIG. 6: (Color online) Critical QS quantities for OCM net- 
works with degree exponent 7 = 2.25, cutoff exponent uj = 
7 — 1, and minimum degree fco = 6. Symbols represent QS 
simulations, solid lines are power regressions and dashed lines 
power laws with exponents ±1/2. 



correction to scaling in g. For 7 = 2.25, the exponents 
obtained with a direct power law regression exhibit an 
excellent agreement with HMF theory, with the excep- 
tion of the scaling of r vs N/g, probably due to finite 
size effects induced by preasymptotic corrections for the 
mean degree (k). 



C. The natural cutoff 

The structural cutoff u = 2 has been used to in- 
vestigate dynamical process H31 [33] in SF networks 
because it prevents degree-degree correlations [25] and 
makes the HMF approach analytically feasible. How- 
ever, the slow growth of k c with N implies large correc- 
tions to the asymptotic scaling [151 I16j and this makes in 
turn the numerical analysis hard in highly heterogeneous 
quenched substrates. For this reason we additionally sim- 
ulate the CP on the ordered configuration model (OCM) 
described in Sec. |ll] with a hard cutoff k c = TV 1 ^ 7-1 ). 
For this cutoff, HMF predicts the critical exponents to be 
v = inax[l/( 7 -l), 1/2] and a = max[(7-2)/( 7 -l), 1/2]. 
We present only the results for highly heterogeneous case 



7 = 2.25 with fco = 6, since the conclusions hold also 
for other values of 7. The crossing points in the mo- 
ment ratio analysis converge faster to the asymptotic 
value, when compared with the UCM. We find a critical 
point p c = 0.4908(1) and moment ratios Mf x = 1.714(4), 
Mfj = 2.29(1), and M% 2 = 3.78(3), which are closer to 
the HMF predictions than those shown in Table |T] for 
UCM networks. As can be seen in Fig. [6j the QS quanti- 
ties exhibit very good scaling properties even for highly 
heterogeneous networks. The critical exponents obtained 
by direct power law regressions are v = 0.774(5) and 
6 = 0.227(5). They differ from the HMF values 0.8 and 
0.2, but also in this case the discrepancy is due to the 
preasymptotic scaling of g. Indeed the regressions includ- 
ing the factor g yield S v = 0.498(4) and S a = 0.502(5), in 
remarkable agreement with HMF. This excellent agree- 
ment is rather surprising since the HMF exponents are 
obtained under the hypothesis of no degree-degree cor- 
relations [HI [TB], while OCM networks are disassorta- 
tive [35]- The present findings suggest that this kind of 
correlations does not strongly affect the behavior of the 
CP. Further theoretical investigation is needed to clar- 
ify why topological correlations play such a minor role in 
this problem. 



IV. CONCLUSIONS 

In this paper we have investigated the behavior of the 
contact process on quenched networks with power-law 
distributed degrees. By performing quasi-stationary sim- 
ulations we have probed the system behavior in the crit- 
ical region, obtaining very precise numerical results for 
both the position of the transition between absorbing and 
active phases and the associated critical exponents. We 
have taken advantage of the recent progress in the ana- 
lytical understanding of the anomalous finite size scaling 
behavior on annealed networks. 

It turns out that the quenched structure of the sub- 
strate has little effect on the contact process dynamics. 
Both the qualitative form of the scaling behavior and 
the quantitative value of the exponents are the same for 
quenched and annealed networks. The only, minor, mod- 
ifications affect the position of the critical point, which is 
shifted in quenched networks due to dynamical correla- 
tions, as well as the values of the moment ratios at crit- 
icality, which depend on 7 on quenched networks, while 
are universal on annealed ones. 

The conclusion that can be drawn is that heteroge- 



neous mean-field theory correctly describes (within our 
error estimates) the critical behavior of the contact pro- 
cess on quenched networks. This conclusion corrects the 
apparent failure of HMF predictions for CP on quenched 
networks that was reported some years ago in Ref. [7] 
With the benefit of hindsight, it is clear now that the 
puzzling discrepancy found in [7] between HMF and nu- 
merical results on quenched networks was not due to sub- 
tle effects of the quenched topology, rather to the non- 
trivial anomalous nature of finite size effects on annealed 
networks. 
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